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Abstract 

According to the delayed-oscillator picture of ENSO, a positive 
SST anomaly in the eastern tropical Pacific will cause westerly wind 
anomalies more close to the dateline that first give a positive feed- 
back, and later, via planetary wave reflection at the western bound- 
ary, a negative feedback. The aim of this study is to follow a chain of 
sensitivities that lead to a delayed-oscillator mechanism in a general 
circulation model. To this end, the adjoint of such an ocean model is 
used for studying sensitivities of ENSO indices. 

The ocean model used in this study is the HOPE ocean general 
circulation model. Its adjoint has been constructed using the Ad- 
joint Model Compiler. Applied to a scalar function computed with a 
forward model run, an adjoint run goes back in time and calculates 
sensitivities as the derivatives of this function to forcing fields or ocean 
state variables at earlier times. 

Results from six adjoint runs are reported, tracing the sensitivities 
of the NIN03 and NIN03.4 indices in October 1987, December 1987 
and December 1988, as simulated by a Pacfic-only version of HOPE 
forced by ECHAM-3 fluxes. 

The sensitivities to sea level can be followed back in time for more 
than a year. They are non-local: patterns propagate back in time that 
are identified as adjoint Kelvin and n = 1, 2 and 3 Rossby waves, with 



speeds compatible with those obtained from model density profiles. 
Both the first and the second baroclinic modes seem to play a role. 
In contrast, the model sensitivities to heat flux, zonal surface currents 
and SST are local and decay in about a month. 

The sensitivities to the wind stress agree with the wave interpre- 
tation of the sea- level sensitivities, but only the n = 1 Rossby wave is 
visible. Going back in time, the sensitivity to westerly anomalies along 
the equator changes sign, in agreement with the delayed-oscillator pic- 
ture. 

Finally, we use a statistical atmosphere model to convert sensitiv- 
ities to wind stress at a given time to sensitivities to SST through the 
atmosphere at that time. Focussing on the sensitivities to the ENSO 
index region itself at an earlier time then closes the circle. These sensi- 
tivities have a natural interpretation as delayed-oscillator coefficients 
and show the expected behaviour of a positive sensitivity in the recent 
past changing to a negative sensitivity at longer lags. However, the 
strength of these feedbacks, and hence the relevance of this mechanism 
in ENSO simulated in HOPE, cannot be determined accurately. 
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1 Introduction 



The first descriptions of the Southern Oscillation were given by pildebrands-| 
son (18971) ; [Walker (1923,1924|) . PJerknes (19"66|) recognized that the periodic 
warming of the eastern tropical Pacific, El Niiio, was part of the same coupled 
system, now called ENSO. He also gave a physical mechanism that identified 
a positive feedback loop that explains the large amplitude of the oscillations. 
The next goal is the identification of the mechanisms behind the changes: the 
onset of a warm (or cold) episode and the negative feedback that limits the 
duration of an episode. This knowledge can be used, if possible, to predict 
ENSO 3-12 months ahead of time. The last decade yielded much improved 
data from the TOGA-TAO experiment and remote sensing, and modelling 
capabilities have increased steadily. 

These models span the range from a stochastically forced harmonic os- 
cillator ( Purgers, 1997|) and canonical correlation analyses ([Barnston and 



Ropelewski, 1992| ) to coupled shallow- water models (e.g., [Zebiak and Cane^ 



19871 ; pchopf and Suarez, 1988] [Kleeman et al., 1995| ) and general circula- 



tion models ( Ji et al., 1996 ; Barnett et al., 1993 ). Several hypotheses for 



the mechanisms that are important in the observed ENSO signal have been 
proposed, and identified in relatively simple ocean-atmosphere models that 
give rise to ENSO-like oscillations but can still be analysed analytically. 

A mechanism that operates in many models is the delayed oscillator 
( Suarez and Schopf, 1988| ; Battisti and Hirst, 1989|) . Consider SST anomalies 
in the eastern equatorial Pacific. On short time scales a positive anomaly 
will weaken the temperature difference with the western Pacific, and hence 
weaken the trade winds, which in turn diminishes the amount of upwelling. 
The subsiding trade winds also generate Kelvin waves that propagate to the 
east and deepen the thermocline there, levelling off the wind-generated slope. 
Both effects lead to even higher temperatures. This classical positive feed- 
back mechanism ( Pjerknes, 1966|) is kept in check by non-linear effects and 
by the effect of an accompanying equatorial Rossby wave, which reflects off 
the western coasts as a negative anomaly Kelvin wave and hence gives a neg- 
ative feedback some time 6 later. The delayed oscillator contains these three 
effects as 

^ = aN-bN^- cN{t - 6) , (1) 

with N some index of the SST in the eastern Pacific, and a,b,c,6 > 0. The 
A^^ term does not influence the period very much, but is the simplest way to 
represent non-linear effects that keep the oscillations bounded. 

This relation has been diagnosed in many models of the tropical Pacific 
( McCreary and Anderson, 1991| ), and gives rise to oscillations with time scales 
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comparable with the observed ENSO record. Evidence that this mechanism 
is present in nature has been presented by [Kessler and McPhaden (1995|) , 
although only the ending of the studied warm event is attributed to western 
boundary reflection. 

However, many other possible modes have been proposed to be active in 
ENSO oscillations. The coupling between ocean and atmosphere may give 
rise to a 'slow SST mode' ( [Neelin, 1991|) , whereas other autors stress the 
importance of the 'recharge' mechanism ( |Wyrtki, 1975| ; |Jin, 1997| ). Also, 



atmospheric events (westerly wind-burst. Madden- Julian oscillations) may 
be responsible for driving ENSO oscillations ( [Kessler et al., 19951 ). A full 
investigation of these mechanisms falls outside the scope of this article, which 
focuses on the role of the ocean. 

In this study we try to identify the delayed oscillator mechanism operating 
in an ocean general circulation model. The approach we take is based on 
the observation that most indicators for ENSO are highly correlated. Many 
effects can be related to a single function, for example, the NIN03 index 
of sea surface temperature in the eastern Pacific. We rephrase the question 
to the cause of ENSO as the computation of the sensitivity of this index to 
changes in the forcing fields or state of the ocean at earlier times. For small 
deviations this question can be addressed with an adjoint model, which can 
trace all influences on one scalar function back in time. 

We describe the HOPE OGCM in section ||, together with a statistical 
atmosphere model that is used in some experiments. This ocean model de- 
scribes ENSO reasonably well (section |2.3|) . In section ^ we describe the 
adjoint ocean model. Using this we search for the sensitivity of a NINO 
index to forcing fields (wind stress, heat flux) and state variables (sea level 
height, zonal currents) in section ^. With the help of the adjoint statistical 
atmosphere model we transform the sensitivities to fluxes into sensitivities to 
SST anomalies, which are then simplified into a description resembling the 
delayed oscillator in section ^. Our conclusions are presented in section |^. 
Two appendices contain technical details and validation runs of the adjoint 
model. 



2 Model Description 
2.1 The ocean model 

The ocean model used in this study is HOPE (The Hamburg Ocean Primi- 
tive Equation Model, |Frey et al., 1997| ; [Wolff et al., 1997| ). As we are mainly 



interested in ENSO we use a domain that is limited to the Pacific Ocean 
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(120°E-70°W, 55°S-67°N). The model has a 2.8° resolution, with the merid- 
ional resolution gradually increased to 0.5° within the region 10°N to 10°S. 
There are 20 vertical levels which are irregularly spaced in such a fashion 
that ten levels are placed in the upper 300 meters. The model contains a 
realistic bottom topography. The numerical scheme uses the Arakawa E-grid 
([Arakawa and Lamb, 1977| ) and a two-hour time step. We use a 360-day year 



subdivided in 12 equal-length months of 4 weeks each. There is no diurnal 
cycle. 

A Newtonian relaxation is employed to restore the surface salinity within 
the whole domain to the climatology of |Levitus (1982[ ) using a time constant 
of 40 days. Otherwise, no explicit fresh water flux is provided. Prognostic 
variables are the three-dimensional fields of horizontal velocity, temperature 
and salinity, and the sea surface elevation. Physical parametrizations in 
HOPE include solar penetration below the surface. About 14% of the so- 
lar radiation incident on the ocean surface is allowed to penetrate beneath 
the 20 m surface layer of the model ([Paulson and Simpson, 1977| ; [Schneider 



et al., 1996| ). The vertical mixing in HOPE is based on a Richardson-number 
dependent formulation and a simple mixed layer scheme to represent the 
effects of wind stirring (see [Latif et al. (199^ ) for details). Further details 
about the model physics as well as the numerical scheme are given in the 
HOPE documentation ( |Woltf et al.,TW70 . 



2.2 The atmosphere model 

For most experiments in this paper we used prescribed heat flux and wind 
stress fields. These were generated with the ECHAM-3 atmosphere model 
from observed SST fields and stored as daily values. ECHAM-3 is the Ham- 
burg version of the European Centre operational weather forecasting model. 
It is described in detail in two reports ( [Roeckner et al., 1992| ; [DKRZ, 1992[) . 



The SSTs were taken from a data set of the British Meteorological Office 
(UKMO), which is referred to as the 'GISST' data set ([Parker et al., 1994[) . 



Additionally, we imposed a relaxation to observed SSTs of 40 Wm~^/K, us- 
ing a Newtonian formulation. This is required to compensate for model drift, 
which produces systematic differences of the simulated SSTs compared to the 
observed ones that were used to generate the atmospheric forcing with the 
ECHAM-3 model. 

A statistical atmosphere model similar to the one used by [Barnett et al 



(1993[) is used to convert the sensitivities of the ocean model to these fluxes 
into sensitivities to SST. This simple atmosphere model is essentially a re- 
gression relationship: wind stress anomalies and surface heat flux anomalies 
are related to sea surface temperature anomalies. This dependence is al- 
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most instantaneous and non-local. It reproduces the large-scale response of 
the atmosphere to large-scale SST anomalies. In particular, the Bjerknes 
positive feedback of reduced trade winds resulting from a reduced east-west 
SST gradient is wcll-rcproduccd in this type of model. On the other hand, 
smaller-scale atmospheric responses are not caught. 

The regression matrix is derived from an integration of the atmospheric 
general circulation model ECHAM-3 with prescribed observed global SSTs 
during the period 1949 to 1990. The linear regression is done in an Empiri- 
cal Orthogonal Function space, retaining the first 10 EOFs of the anomalies, 
which account for 63% and 49% of the field variance in the SST and the at- 
mospheric fields, respectively. The dominant modes of interannual variability 
are found to be well captured by the EOF truncation. To compensate for 
systematic deficiencies of the OGCM, the AGCM's simulated SST anomaly 
is corrected using another regression matrix based on the leading 10 EOFs 
of observed and simulated SST anomaly. 

Let the observed and simulated SST anomalies be denoted by Tohs{x,t) 
and Tsira{x,t), respectively. The zonal and meridional wind stress anomalies 
were pooled together with the surface heat fiiix anomalies in one data set of 
atmospheric fiux anomalies denoted by ^echaml^:^, The EOF expansions of 
the observed and simulated fields are 

Tohs{x,t) = ^Q;„(t)e„(x) (2) 

n 

m 

TsUx,t) = Y.'yiit)9i{x) . (4) 

I 

The linear regressions yield matrices C-^-* and C*-^-* whose entries are given 

by 

Cff = (6) 

Angle brackets signify time averages. The estimated wind stress anomaly, 
J^esti^, t), may now be obtained by the simulated SST anomaly via 

m 

with 

Pe.t,m{t)=j:C!:iC^hlit)- (8) 
n.l 
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SST anomalies simulated by the HOPE model when forced with the 
ECHAM-3 data mentioned above are not identical to those observed (Fig. |]). 
To compensate for these systematic differences, the SST anomalies produced 
by the HOPE model are corrected before they are used to infer the atmo- 
spheric fields. For this correction another regression matrix based on the 
leading 10 EOFs of simulated and observed SST anomalies is used. A simi- 
lar approach was taken by p3arnett et al. (1993|) . The additional surface heat 
flux correction is given by 

Qcorix,t) = \{T^imix,t) -Tohs,estix,t)) (9) 

Tobs,est(x,t) = EE^i^^Ki)e„(x), (10) 

n I 

with the relaxation constant A = — 10 Wm~^/K. As the estimates for the 
wind stress and surface heat flux anomalies are always made for the next 
month to be integrated, the linear regressions are computed with the SSTs 
leading the atmospheric fields by one month. 

In Fig. H we compare the zonal wind stress along the equator gener- 
ated with the statistical atmosphere model with the flux generated by the 
ECHAM-3 model, which was used as a forcing field. The large-scale features 
are seen to correspond quite well. 



2.3 ENSO simulated by HOPE 

Next we establish whether the model used is in fact able to reproduce the 
observed behaviour of ENSO. The ocean model was spun up from the Levitus 
climatology for 30 years with climatological forcing. It was then integrated for 
42 years starting in January 1949, using wind stress and surface heat fluxes 



from the ECHAM-3 atmosphere model integration described in section ^ 
as upper boundary conditions, including the 40Wm~^ relaxation to observed 
SST. 

The spatial structure as well as the temporal evolution of observed and 
simulated interannual SST anomalies in the Pacific Ocean can be inferred 
from the pattern of the first EOF and timeseries of the associated principal 
component. Overall, the anomaly structure, shown in Fig. |^, is simulated 
reasonably well by the HOPE model. There are, however, some differences 
to the observations. The simulated SST anomalies are too much equatori- 
ally confined, and they extend too far into the western Pacific. The latter 
deficiency is related to the too strong cold tongue which extends too far into 
the western Pacific. Hence, in addition to the average SST in the NIN03 
region (5°S-5°N, 90°W-150°W) we will use a slightly further western region 
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Figure 1: The zonal wind stress in Nm~^ averaged from 5°S to 5°N gen- 
erated from observed SST by the ECHAM-3 model (a) and from simulated 
SST by the statistical atmosphere model (b). 
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Figure 2: Maps of the first EOF of monthly mean SST anomahes from 1949 
to 1990 as derived from (a) the HOPE model simulation and (b) observations 
(GISST data set) explaining 26% and 29% of the total variance, respectively. 
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Figure 3: Time series of the principal components. 



(NIN03.4: 5°S-5°N, 120°W-170°W) as an ENSO index in our sensitivity 
studies. The comparison of the principal components, shown in Fig. ^, re- 
veals satisfactory agreement (r = 0.88) between the temporal evolution of 
the simulated and observed SST anomalies. 

The HOPE model can be coupled to the statistical atmosphere described 
in the previous section instead of being forced by prescribed upper boundary 
conditions. Such a combination of an OGCM and a statistical atmosphere 
model is generally referred to as a Hybrid Coupled Model (HCM). It shows 
a regular oscillation of the SST anomalies with a period of 2.5 years. The 
spatial structure and amplitude of the SST anomalies are comparable to 
those obtained with prescribed forcing (Fig. ^ja). 



3 The adjoint model 

In this section we define what we mean by the term 'adjoint model' and 
present some aspects of the adjoint of the OGCM we constructed. 
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3.1 Definition 



The system being studied is described by a state vector x of prognostic 
variables. Let the model A4 be some function, in general non-linear, that 
computes from an initial state xo and a set of forcing fields yi at times 
ti {i — 1, . . . ,n) a state vector Xn at time tn 

Xn = M{xo,yi) . (11) 

The tangent model M{xo,yi) is the linear mapping that gives the effect 
on Xn of infinitesimal changes in the initial state Sxq and in the set of forcing 
fields Syi 

X A^i \ r M(xo + e5xo,yi + e5yi) - M{xo,yi) 
dxn^ M{xo,yi) j =hm . (12) 

Let f{xn) be some scalar function of the final state, for instance an averaged 
field, index or cost function. We can now compute the change in this function 
of the final state as a result of a change in the initial state and forcing fields 
as 

with the mapping M^{xo,yi) the adjoint model. It propagates the sensitivi- 
ties of the scalar / to the final state back to the initial state, and can be seen 
as a backward derivative of the model. 

In practice, the model M. consists of a series of discrete time steps Xi — > 
Xi+i{xi,yi) for i = 0, . . . ,n - 1, so 

•^n — •^n('^n— 1) Z/n— l) 

— ^n(^n-l(^n-2, yn-2), Z/n-l) 

(14) 

The adjoint model is then evaluated by time-stepping the discrete model 
back from t„ to to, at each step applying the chain rule of differentiation 
with respect to the model variables 

dfjxn) ^ df{xn) dXj+i 
dxi dxi+i dxi 
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and the forcing fields 

dyi dxi+i dyi 

Note that the evolution of the adjoint fields Of /dxi does not depend on the 
function /, hence except for the initialization of df /dxn the adjoint model 
is independent of the function /. Physically, the adjoint field represents the 
sensitivity of / to infinitesimal perturbations of the prior state Xi and forcing 
field yi. 

For a non-linear model M., the response to a finite perturbation may or 
may not be similar. An example of the latter case is convection, which is 
often represented by a step function of the difference in density of two layers. 
The derivative of this function to the density is zero almost everywhere, yet a 
finite perturbation can induce or inhibit convection. The effect of non-linear 
terms in a particular problem can be gauged to some extent by comparing 
the effect of a finite perurbation 5xo, 5yi on the scalar function f{xn) to the 
effect deduced from the derivatives, Eqs (|12], pT^). 

The two types of adjoint fields in Eqs (|T5|) and (|T6|) are conceptually dif- 
ferent. Derivatives to forcing fields df /dyi (for instance externally prescribed 
fluxes) are not used at times t < ti. These adjoint fields thus represent quan- 
tities that can be varied independently. On the other hand, derivatives to 
internal model state fields df/dxi are used in the previous time step to com- 
pute df/dxi-i, and hence cannot explicitly be varied without invalidating 
the adjoint model for t < ti. 

As defined above the adjoint field is an extensive variable: the influence of 
a single grid box is proportional to its size. Also, the influence is proportional 
to the time step used. We therefore normalize the adjoint field to the grid 
box size and a time interval of one day and give the derivative to forcing 
fields (fluxes) in units [/]/([a;] day sr). 

3.2 Adjoint HOPE and atmosphere models 

We constructed adjoint models (as defined above) of the ocean and atmo- 
sphere models in the following way. 

In HOPE, a time step consists of a cycle of separate tasks. For each 
of these subroutines (or easily definable subdivisions) the adjoint code was 
generated using the Adjoint Model Compiler ( |Giering and Kaminski, 1996| ; 
Giering, 1994| ). This tool gives computer code which is the adjoint of the orig- 
inal algorithm. These adjoint routines were hand-optimized in some cases 
and integrated into an adjoint model. Further details can be found in Ap- 
pendix 0. 
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This construction gives the exact adjoint of the complete (discrete) model. 
We realize that there are many partial derivatives that can probably be safely 
neglected, like the dependence of the mixing parameters on the temperature 
and the detailed intricacies of the equation of state. However, it is a nontrivial 
task to verify that omissions are in fact harmless, and we did not (yet) 
perform this analysis. 

There is one point where it was necessary to deviate from the exact deriva- 
tives to include some non-linear effects. In our model the vertical eddy vis- 
cosity Ay and diffusivity Dy depend on the Richardson number R as 

= -^Y^^^ + other terms (17) 

Dv = (1^;^ + other terms (18) 

R = -S^P/^' (19) 

R! = max{0,R) (20) 

with Ca = Cd = 5. For marginally stable layers (0 < i? < i?cut = 0.3) 
the large derivatives dAy/dR, dDy/dR are not very good estimates for the 
effect of a finite perturbation, as the functions have very sharp kinks at 
= 0. In practice this manifests itself as dipoles of very strong positive and 
negative influence one vertical grid point apart in the adjoint temperature 
and salinity flelds. This instability can be suppressed by smoothing the 
derivatives of Eqs (p!7|-[l8D. In the evaluation of the derivative, we therefore 
use R' = RI^JR for < i? < i?cut- 

In contrast, the non-linear terms associated with convection were found 
not to influence the validity of adjoint runs with targets the NIN03 and 
NIN03.4 indices. This was verifled by using different-sign perturbations in 
the flnite difference experiments of App. 0, which show effects of less than 
35%. The physical reason is that the signal region, in the cold tongue, is 
usually stably stratifled, so that convection does not play a role there. Signals 
to the target region often move in the thermocline, where convection also is 
irrelevant. 

A recent article ( [Sirkes and 'Iziperman, 1997| ) cautions against using the 



adjoint of the finite differences model to trace sensitivities in time, as a 
discretization that is stable in the forward runs may be unstable in the ad- 
joint model. The specific example shown, a leap-frog scheme, is not used in 
HOPE, and we did not find this kind of instabilities in our adjoint runs in 
the variables shown. There is a slight two-point instability in ON^/dTy in the 
equatorial western Pacific that may be caused by similar mechanisms. 
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The adjoint of the statistical atmosphere model was constructed in a 
similar fashion. As this is a linear model, it does not depend on the forward 
fields, which simplifies the code. The time lag between the SST and fluxes 
that the atmosphere model introduces has been taken into account. 



4 Sensitivity experiments 

To investigate the causes of temperature changes in the eastern tropical Pa- 
cific in the HOPE model we performed a few experiments with the adjoint 
HOPE model described in section |3]^ using prescribed fluxes. These experi- 
ments are the exact opposite of perturbation experiments that search for all 
effects of one given perturbation, usually a wind burst (e.g., piese and Harri 



son, 1990| ). We search backwards in time to find all sensitivities of one target 



function, which we take to be the NIN03 (or NIN03.4) index of SST in the 
eastern Pacific. A sensitivity field gives the change in the target function if 
we would have added a perturbation at that time and place. 

From perturbation studies and lag-correlation analyses (|Gill, 1982| ; phi- 



lander, 1990| ) it has been found that many perturbations propagate as equa- 



torial Kelvin or Rossby waves. The baroclinic Kelvin waves travel eastwards 
within a few degrees of the equator, with a speed of 2-3 ms~^. The n-th mode 
Rossby waves travel to the west at l/{2n + 1) times this speed, with a maxi- 
mum perturbation to the sea level for the first three modes at about 5°, 8-10° 



and 12-15° off the equator (in theory, these are y3/2, and y (5 -|- -\/l8)/2 
Rossby radii). This propagation implies that sensitivities propagate back- 
wards along the same trajectories, which we refer to as adjoint Kelvin or 
Rossby waves. Reflections of these waves at the coast are therefore the back- 
ward analogon of reflection of the usual waves: at the western boundary, an 
adjoint Kelvin wave coming from the west will reflect into an opposite sign 
adjoint Rossby wave coming from the east. We will verify this interpretation 
by checking that the propagation speed is consistent with predictions from 
linear theory. 

As background for the sensitivity experiments we used the years 1987 
and 1988. In 1987 there was an El-Nino episode (the observed NIN03 index 
was -1-1.5 K to -1-2.0 K for much of the year, but had returned to 1.0 K 
by the end of December), whereas 1988 ended on a strong cold anomaly 
(—1.5 K). As target dates we took October 1987 (warm), December 1987 
(transition) and December 1988 (cold). These three dates give an impression 
of the dependance of the sensitivities to the background state. For a linear 
model the three runs should give identical results. Our resources were not 
sufficient to do a systematic study of these non-linear effects. We studied 
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the sensitivity of the NIN03 index N3 (NIN03.4 index A^3,4) at these times 
to various fields at earher times. The interesting ones are the input fiuxes 
(wind stress r^^y and the heat fiux Q) and state variable sea level height (. 
The former have been averaged over 7.5 day 'weeks', the latter is sampled 
at the beginning of these periods^. Unless otherwise noticed the illustrations 
are from the run to the NIN03 index at the end of December 1987. 



4.1 Sea level 

The clearest signals can be seen in the sensitivity to the sea level. This is due 
to the fact that, as in the forward case, the adjoint Kelvin and Rossby waves 
are clearly separated by latitude. The interpretation of the adjoint field is 
however complicated, as there is no physical way to change the value of this 
state variable directly. A balanced change in the sea level requires changes 
to the subsurface structure as well, and the effect of these are not included 
in the partial derivatives shown here. 

The sensitivity of the SST to the sea level on short time scales is caused 
by the following chain of infiuences. An increase in sea level will depress the 
thermocline slightly^. This deeper thermocline in an upwelling region in turn 
causes the sea surface to be less cold. It is therefore no surprise to see that 
at the end of December 1987 the NIN03 index would have been higher if 
the sea level would have been higher one week earlier in the same area. The 
NIN03 index is the average temperature in the middle box shown in Fig. ^a. 
The sensitivity is thus positive and confined to this region, with an emphasis 
on the upwelling region. (The slightly negative field in the western Pacific is 
a weak barotropic mode.) 

Going back in time, at the beginning of November (Fig. Qd), one can see 
an adjoint Kelvin wave coming from the western Pacific, and some adjoint 
n = 1 Rossby waves coming from the eastern boundary. The adjoint Kelvin 
wave should be interpreted as follows: if the sea level would have been (in- 
finitesimally) higher in the vertically hatched area in November 1987, then 
there would have been a baroclinic Kelvin wave added to the actual state 
at the position indicated, and the temperature would have been (infinites- 
imally) higher in the NIN03 area 2 months later. The Kelvin wave would 
have depressed the thermocline in that region. (Note that this wave is not 
actually added to the state: the linearization is performed around the actual 
background, not around climatology.) 



^Mpeg movies and space-time diagrams of the experiments are available at 



tittp : //www.knmi .nl/^oldenbor /Sense/pictures .html. 

■^Most extra water will just flow away though (Fischer and Latif, 1995) 
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Figure 4: The sensitivity of the NIN03 index to sea level changes dN^/dC, in 
K/m/sr at the end of 1987 in the last week of December (a), the first weeks 
of November (b), July (c), April (d) and January 1987 (e), and July 1986 
(f). The boxes denote, from east to west, the NIN012, NIN03 and NIN04 
index regions. -^^ 



In the beginning of July (Fig. all the adjoint Rossby waves going to 
the measurement region from the eastern boundary have reached the coast. 
Comparing with the adjoint run of a linear shallow water model (not shown), 
it seems that these are partially reflected into an adjoint Kelvin wave, which 
merges with the tail of the original one. The resolution does not permit 
coastal Kelvin waves to travel north and south. On the other side, the 
adjoint equatorial Kelvin wave has reflected at the western coasts into a set of 
adjoint Rossby waves of opposite sign. The first blob at 165°W corresponds 
to the n = 1 reflection off New Guinea (cf. Fig. |), the second one (at 
150°E) has been reflected off the Philippines. (The HOPE land-sea mask 
has approximately two north-south walls at these locations.) At the same 
longitude the n = 2 reflection off New Guinea is forming at 8°S, and further 
west the beginnings of higher-order adjoint Rossby waves can be discerned. 
The presence of the parity-odd n = 2 mode is due to the fact that the western 
coast has some structure in our model. 

Three months earlier (Fig. ^) the first adjoint Rossby wave has almost 
dissipated, the strongest signal is now at 160°W. The southern n = 2 signal 
is now at the dateline, whereas the northern reflection has reached 150°E at 
8°N, clearly separated from the n = 3 reflection at 15°-20°N. The situation 
at the beginning of January is not much different (Fig. ^). 

At the beginning of July 1986, 18 months before the NIN03 index is 
evaluated (Fig. ^), the only influence remaining is a pattern of higher order 
Rossby waves that slowly traverses the basin. This sensitivity on the sea level, 
and equivalently the thermocline depth, of the off-equatorial warm pool is in 
agreement with the original recharge hypothesis ( |Wyrtki, 1975|) . However, 
a sensitivity of the NINO indices to the zonally-averaged thermocline depth 
along the equator at this time scale ( Jin, 1997 ) is not seen. 

The reflection pattern off the western coasts is clearly visible in the space- 
time diagrams of 5° latitude bands at 5°N, the equator and 5°S in Fig. |^. In 
1987 the adjoint Kelvin wave, travelling at about 2ms~^ (solid line in Fig. |^) 
reaches the coast of New Guinea in the first week of October, generating a 
Rossby wave at 5°S that picks up a northern component some time later. The 
main reflection at 120°E occurs at the end of September. In 1988 the adjoint 
Kelvin wave is both faster and weaker, arriving at the eastern boundary 
one week earlier. East of 150°W there are also the adjoint Rossby waves 
generated by the edges of the NIN03 region, and adjoint Rossby waves with 
opposite sign caused by the adjoint Kelvin wave leaving the NIN03 signal 
region. They reach the eastern coast within a few months. These waves are 
much stronger in the cold conditions of 1988. In HOPE the NIN03 region is 
characterized by intense currents and changes in thermocline depth at this 
time, so that linear theory is not applicable. 
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Figure 5: The sensitivity of the NIN03 index to sea level changes in December 
1987 dN^/dC, in K/m/sr at 5°N (a), the equator (b) and 5°S (c), averaged 
over 5° latitude bands; the same one year later, to December 1988 (d,e,f). 
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A second adjoint Kelvin wave is present in the western and central Pacific 
in August and July. It seems to be generated from three sources: partial 
reflections at the eastern boundary (especially in the 1988 experiment with 
strong activity in the eastern Pacific), a partial reflection of the first adjoint 
Rossby wave at the edge of the warm pool (mainly in 1987), and finally a 
reinforcement near the western boundary by a reflection off the western coast 
of New Guinea by the Rossby wave reflected off the Philippines. The arrival 
of this second Kelvin wave generates a third adjoint Rossby wave during July. 
Similar effects seem to produce two more adjoint Kelvin and Rossby waves 
in February and November the year before respectively. 

These reflected Rossby waves merge into a structure with an effective 
group velocity that is significantly lower than the phase velocity of the indi- 
vidual n = 1 Rossby waves, indicated by the dashed lines in Fig. ^,c with 
speed 0.6 ms~^. This behaviour is very similar to the properties of the ob- 
served Kelvin and Rossby waves that constitute ENSO oscillations ( |(Jhao 
and Philander, 1993| ; [Kessler et al., 1995|) . 

These reflection patterns are compatible with observations from sea level 
measurements near the eastern and western boundaries (Poulanger and Fu. 
1996| ). Some of the discrepancies at the western coast reported in this paper 
may be explained by the multiple reflections observed in our model. 

4.2 Wave speeds 

To check the interpretation of Kelvin and Rossby waves, we have measured 
the speeds of the structures visible in the adjoint run and compared the 
results to wave speeds calculated from the density structure of the model. 

The speed was estimated by taking 1- and 2-month lag correlations for 
each latitude. From Fig. |^ one can read off the following speeds for the 
waves. In the last month of the run leading to the NIN03 index at the end 
of December 1987, the equatorial Kelvin wave has a speed that varies between 
1.8 ms'^ at rS, 2.0 ms"^ on the equator and 2.4 ms"^ at 2°N (Fig. |a). In 
the experiment leading to October 1987 we find speeds of 2.1 to 2.9 ms^^ 
(Fig. Pd), and at the end of 1988 they are around 3 ms~^ (Fig- l^c). 

The Rossby waves during these times in the Eastern Pacific have speeds 
at 5°N of 0.5 ms~^, 0.5 ms~^ and 0.7 ms~^, this is somewhat slower than 1/3 
the speed of the Kelvin waves due to the shallower thermocline. The speeds 
of reflected adjoint Rossby waves in the central Pacific have been read off 
Fig. pi. The northern branches of the Rossby waves move with —0.55 ms~^ 
(5°N), -0.34 ms"^ (10°N) and -0.15 ms"^ (19°N). The southern lobes are 
slower and further from the equator: —0.45 ms~^ (5°S), —0.23 ms~^ at 12°S 
and — 0.12 ms"^ (19°S, only visible in 1986). These speeds are roughly the 
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Figure 6: The 1-month lag correlations p as a function of latitude and speed 
in m/s of dN^/dC, in the last month of the Dec 1987 run (a), the Oct 1987 run 
(b), the Dec 1988 run (c), and the 2-month lag correlation in January-April 
1987 (d). The contours refer to — log(l — p), the black dots and hnes refer 
to the values mentioned in the text. 
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same in the other two experiments, except that the n = 1 mode is often 
poorly visible. 

The observed adjoint Kelvin waves therefore seem to be a mixture of 
the first baroclinic mode, with speeds around 3 ms~^, and the second baro- 
clinic mode, with a speed of around 2 ms~^. This is in agreement with 
other theoretical studies (e.g., piese and Harrison, 19901) . Observations of 



Kelvin waves give a single mode with a speed of 2.4 ± 0.3 ms ^ from TOGA- 
TAO subsurface temperatures ( [Kessler et al., 1995| ) and 2.9 ± 0.9 ms"^ from 



TOPEX/POSEIDON altimetry data (|Boulanger and Menkes, T995| ; [Chelton 



and Schlax, 1996| ). Also, in these observations the shallowing of the thermo- 
cline towards the east is compensated by a stronger density contrast across 
the thermocline, leading to an almost constant speed up to 110°W. In con- 
trast, in HOPE we observe lower speeds east of 160°W (see Fig. as the 
thermocline is not sharp enough due to discretization effects. 

These speeds can also be computed from the density profile^ (see, e.g.. 
Gill, 1982| ; [Philander, 1990| ), using a rigid lid to eliminate the barotropic 



mode. The first three modes of the horizontal velocity are shown in Fig. |^. 
The first baroclinic mode has a zero around 1500 m and a Kelvin wave speed 
of 2.6-3.3 ms~^, in agreement with the profiles discussed in [Philander (1990D ; 



piese and Harrison (1990] ). This mode seems to be associated with the rise in 
salinity around that depth in the Pacific, which is absent in the Atlantic. The 
next modes have their first zeroes near the thermocline, at around 100-150 
m in the equatorial zone, substantially deeper elsewhere. The speed fields 
are shown in Fig. R|. These maps can be compared with the ones obtained 



from much higher-resolution data in picaut and Sombarddier (1993[ ). They 



find speeds that are 20-25% lower than the ones shown here, whereas [Kessler[ 



and McPhaden (1995[) use a speed measured at one point (150°-158°W) that 
is about 10% lower than the value we find there. We assume the differences 
are due to the coarse vertical grid we employed, which is just the 20-level 
HOPE grid. 

With a well-developed thermocline, as in the western and central Pacific, 
and in the east during the December 1987 run (Fig. |^,b), a Kelvin wave in 
the first baroclinic mode mode does not influence the NIN03 index of surface 
temperatures as much as one in the second mode. The second mode changes 
the thermocline depth, which most effectively affects the surface temperature 
in the upwelling regions (see section [IT]] ). Also, the eigenfunction of the first 
mode has a smaller amplitude at the surface so that a zonal stress at the sur- 



■^This computation depends somewhat on the mixing scheme, and is done for mixing 
parameters proportional to which differs from the Richardson- number based scheme 

with a mixed layer of the HOPE model. The difference will be mainly in the absence of a 
mixed layer in the computed modes. 



19 





Figure 7: Profiles of tlie first three baroclinic modes and temperature at the 
equator, 180° (a) and 100°W (b) in Dec 1987, and at 100°W in Oct 1987 (c) 
and Dec 1988 (d) in the HOPE model 
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face does not excite this mode as much as the second one in the derivatives to 
wind stress (section [4.4|) . Finally, the Richardson-number dependent mixing 
scheme, Eqs (|T^-|TB|), tends to decouple the ocean below the thermocline, 
also preferring the second mode. 

However, when the thermocline is absent in the cold tongue, as our model 
generates in December 1988 (Fig. |^), the first mode is most effective in 
changing the SST in the NIN03 area, although the adjoint Kelvin wave is 
weaker in this case (see Figs ^,e). The October 1987 run, with a weak 
thermocline, excites adjoint Kelvin waves in both modes. 

The speed of the second mode Kelvin wave in the December 1987 run 
(Fig. ^) is slightly higher than the 1.8 to 2.1 ms^^ in the area west of the 
NIN03 box in Fig. In 1988 the speed corresponds well with the little 
over 3 ms~^ read off from Fig. j^a. The n = 1 reflected adjoint Rossby waves 
correspond reasonably well with 1/3 the zonally averaged speed of the second 
mode (0.57 ms~^ at 5°S, 0.62 ms~^ at 5°N); the same holds for n = 2 (1/5) 
at 10°N (0.34 ms~^) in all experiments. Note that this means that during 
the reflection off the western boundary energy can flow from the second into 
the first baroclinic mode; the WKB approximation breaks down when the 
equivalent depth changes significantly on scales comparable to the Rossby 
radius. 

The speeds at the other extrema in Fig. |^, corresponding to higher-order 
Rossby modes, are much lower and seem to resemble more the speed of the 
third baroclinic mode (0.21 ms~^ at 11°S, 0.15ms~^ at 19°S and 0.17ms~^ at 
15°N). Unfortunately the derivatives to the horizontal speeds are too noisy 
to determine their vertical structure, we can only assume that the second 
and third modes can also be mixed up in the shallow water at the western 
coasts. 



4.3 Zonal currents 

According to some publications (|Picaut and Delcroix, 199"5| ; [Picaut et al. 



199(j|) zonal currents converging at the edge of the warm pool play a major 



role in ENSO. As can be seen in Fig. ^ this mechanism is not visible in 
the adjoint HOPE model: the NIN03 index only depends significantly on 
surface currents in the measurement region over the last half year. These 
sensitivities seem to correspond to the surface currents of the adjoint Rossby 
waves visible in Figs ^,c in this area. The shallow thermocline in the cold 
tongue may increase the effect of changing just the surface current. 
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Figure 8: The speeds in ms ^ of the first (a) and second (b) barochnic modes 
deduced from the density profiles at the end of 1987. 
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Figure 9: The sensitivity of the NIN03 index to zonal surface currents 
dNs/duo in K/ms~^/sr averaged over 4°S to 4°N. 
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Figure 10: The sensitivity of the NIN03 index at the end of 1987 to zonal 
wind stress, ON^/dTx in K/Nm~^/sr/day, at the end of 1987 (a), in the first 
weeks of November (b), July (c), April (d) and January 1987 (e), and July 
1986 (f) 
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4.4 Wind stress 



The derivative to the wind stress is a superposition of the wind fields as- 
sociated with the Kelvin and Rossby waves. In Fig. we show the situ- 
ation at the same times as the sensitivity to the sea level height in Fig. ^. 
These fields have the following interpretation: if there would have been a 
positive (eastward) anomalous wind field at the position of a positive deriva- 
tive, then the NIN03 index would have been higher at the end of 1987. 
To give a quantitative example: suppose there would have been an extra 
westerly windstress of 0.05 Nm-^ at 180°-160°W, 2.8°S-2.8°N (smoothed 
at the boundaries) during the first eight days of November. The sensitiv- 
ity to the zonal windstress in this region is about 6 K/Nm~^/sr/day, so 
the effect would be an increase of the NIN03 index two months later of 
0.05 Nm-2 X 8days x 20°5.6°7r7(180°)2sr x 6K/Nm-Vsr/day = 0.08 K. Us- 
ing the exact dotproduct between the sensitivity and the extra stress the 
prediction is 0.064 K, a perturbation experiment gave 0.077 K. Similarly, in 
April the sensitivity in this region is roughly —2 K/Nm~^/sr/day, so that 
the same extra wind stress would have given an increase in the NIN03 index 
of —0.03 K nine months later. (Appendix |B| gives comparisons between the 
linear sensitivities and full perturbation experiments.) 

During the first week of the adjoint run (Fig. |TO|a), one sees that a positive 
wind anomaly (weaker trade winds) reduces upwelling, and hence heats up 
the NIN03 region. At the extreme eastern end of the region a positive 
anomaly causes a more divergent wind field, and hence more upwelling and 
a lower NIN03 index. 

At earlier times the disturbance caused by an anomalous wind field has to 
propagate to the eastern Pacific in order to infiuence this index. Comparing 
with the derivatives to the sea level (Fig. ^), one sees the textbook relations 
(e.g., |Gill, 1982D that for a Kelvin wave the zonal wind field and sea level 



have the same shape, whereas an n = 1 Rossby wave has an anti-cyclonal 
circulation around a positive sea level anomaly. The n = 2 wave clearly 
develops its antisymmetric wind field (for instance at 160°E in Fig. p!0|c), 
whereas the n = 3 Rossby wave is symmetric again. However, the wind-field 
sensitivity of the higher order modes is much weaker than the sensitivity 
of the n = 1 mode. The 1988 experiment, during cold conditions, shows 
a weaker dependence on wind stress; discussion of this effect is delayed to 
section ^. 

The result of these adjoint Kelvin and Rossby waves is that during the 
half year leading up to the measurement the sensitivity to at the date line 
is positive, meaning that a weaker easterly trade wind or even westerly wind 
will increase the SST at the NIN03 region. At earlier times, the sensitivity 
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has the opposite sign — a westerly anomaly will cool the NIN03 region at 
these lead times. 

The sensitivity to the meridional wind stress Ty is much smaller, and 
mainly shows the edges of the structures visible in the derivative to the sea 
level, Fig. |. 

4.5 Heat flux 




JASON 



Figure 11: The sensitivity of the NIN03 index to the heat flux dN^/dQ in 
K/Wm~^/sr for various longitudes in the eastern Pacific, averaged over 5°S 
to 5°N. 

Finally, an obvious influence on the surface temperature in the NIN03 
region is the heat flux into this area, or, almost equivalently, the temperature 
some time previously due to persistence. This is relevant only for the last 



few months preceding the measurement. In Fig. |TT]we plot this derivative as 
a function of time at the equator. Except for some waves that are excited in 
August-September the sensitivity drops off exponentially with a time scale 
of about one month (25 ± 5 days). There is some advection from the north- 
east, but this does not extend beyond a few degrees of latitude outside the 
NIN03 area. 

n = 1 Rossby waves are excited at the point where the thermocline is very 
shallow in the HOPE model (around 140°W). We assume that heat input at 
this point lowers the thermocline and thus excites waves. Only the Rossby 
waves travel in the correct direction to influence the NIN03 index. In other 
circumstances a heat flux does not efficiently excite long waves. 
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5 The delayed oscillator in HOPE 



The experiments described in the previous section give the sensitivities of 
an uncoupled ocean model driven by prescribed fluxes. We would like to 
compare this with the concept of the delayed oscillator, which is a coupled 
mechanism. However, the delayed action term is due to long waves prop- 
agating in the ocean, and is already visible in the sensitivity to the zonal 
wind stress changing sign at the equator in Fig. |10[ The necessary feature 
of an accompanying atmosphere model is a large-scale and quick adjustment 
of Tx. to zonal SST gradients along the equator. This can be represented by 
the simple statistical atmosphere model described in |^.2| , which generates a 
zonal wind stress similar to the forcing flelds (see Fig. |ID. 

Using the adjoint of the statistical atmosphere model mentioned in section 



3l2| we can relate the sensitivities of the NIN03 or NIN03.4 index to the 



forcing flelds Ti (the adjoint flelds dN^jri and dNn/dQi we presented in 



sections 0| and |475|) to sensitivities to SST 



months i ^-^^ i OJ-iOli^i 



The flelds dN^/dTi dTi/dTi _i are plotted in Fig. |l^ for the months corre- 
sponding to Fig. |10| (note the linear scale). These are again a kind of forcing 
flelds: as we did not run in coupled mode they are not related to the actual 
SST flelds and only serve to change the fluxes. Hence the summation over 
previous months in Eq. (|2T|) . 

The largest influences are along the Peruvian coast and off southern Mex- 
ico. The latter is a consequence of a mismatch between the version of HOPE 
we used and the version to which the statistical atmosphere was tuned, and 
should be discarded. The runs with as target the NIN03.4 index should be 
less sensitive to this effect. We also flnd large influences of the extra-tropics, 
probably due to spurious correlations in the statistical atmosphere model. 
The sensitivity to ocean temperatures in the cold tongue gives the effect 
of the weakening and strengthening of the Walker circulation. Up to May 



(Figs |T^-c) the influence of heating the cold tongue is positive: the warmer 
water suppresses the trade winds, exciting downwelling Kelvin waves. At 
longer lead times the effect is negative (Figs 0d,e): the same change in 
trade winds also gives rise to downwelling Rossby waves that reflect at the 
western boundary as a series of upwelling Kelwin waves. During the last 
month (Fig. |l2|a) there also is a local negative damping feedback, and the 
effect of warmer waters west of the NIN03 region is to enhance the trade 
winds and decrease the index. 
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Figure 12: The sensitivity of the NIN03 index for forcing the statistical 
atmosphere with an extra SST field dN^/ dTidTi/ dT in K/K/sr in November 
1987 (a), October (b), June (c), March (d) and December 1986 (e) 
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This picture can be reduced to a one-parameter dependence of the NIN03 
index on itself at earher times. Temperature changes in the eastern Pacific 
are highly correlated, and also the sensitivities to these changes in the eastern 
Pacific are due to the ENSO pattern. We expand Eq. (|2T|) in EOFs and 
keep only the first one Ci (shown in Fig. §a), 

= E E e. ■ 5T^., - E -1 • 5T,_, . (22) 

This truncation underestimates the variations by a factor C ~ 1.5. Next we 
realize that all indicators of sea temperature in the eastern Pacific are highly 
correlated, and approximate 

ei ■ 5Ti ^ C(^) 5N, (23) 

■ 7^ « C-^' E J. . (24) 



' index ' 



region 



Over the 42- year simulation period discussed in section ^]3| we find that 5Ti-ei 
has a correlation of r = 0.95 with the NIN03 index, and r = 0.96 with the 
NIN03.4 index; the corresponding constants are C(^^ = 0.57 and = 0.51. 
The derivatives were fitted over the last 12 months of the six experiments. 
Except for the final two months, when the adjoint fields are still close to 
the index region, the agreement is very good, as the adjoint fields are mainly 
confined to the tropical region. The proportionality constants are C^^ = 2.8, 
= 4.2. 

This finally gives the dependence of the index to itself some time earlier 

ON. 



6N^ ^ CC^'^C^'^ Y: ^SN, . (25) 

i=l 



The factor CC(i)C(2) is about 2-3 for the NIN03 index, 3-4 for the NIN03.4 
index. The sensitivities dNn/dNi are shown in Fig. [T3 for the six experiments 



we have performed (to the NIN03 and NIN03.4 indices in October 87, De- 
cember 87 and December 88). The interpretation of Eq. ( p5|) is not entirely 
straightforward, as a change SNi not only causes a change 6Nn directly, but 
also through intermediate steps SNj with i < j < n, which are not explicitly 
mentioned in this form. This way the coupled response to a perturbation at 
6 months' lead time will not only consist of the negative feedback indicated 
in Fig. |T^, but also of the square of the 3-month positive response, and con- 
tributions from all other partitions of six. The effect is to move the negative 
feedback to longer times, in better agreement with the NIN03, NIN03.4 
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Figure 13: The sensitivity of the NIN03 and NIN03.4 indices to the same 
index at the end of 1987 and 1988 through the statistical atmosphere in 
K/K/month. Including the sensitivity to correlated SST changes would 
roughly result in a factor three increase. 
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autocorrelation functions. In principle we could have computed this coupled 
response function from Fig. |13|, but the large uncertainty in the scale make 
this a futile exercise. 

On the other hand, Fig. |1^ and Eq. (p5D can be compared directly to a 
discretized form of the delayed oscillator equation Eq. (P 

Nn = {l + aAt)Nn-i-bAtN^_,-cAtNn-5 , (26) 
<W„ = ((1 + a At) - 3bAtN^_^)6Nn^i - cAt dN^s , (27) 

with At one month. As expected from this equation, the 'immediate' positive 
feedback differs very much between a situation where the index is falling to 
zero (Dec 1987) and a situation in which the index is strongly positive (Oct 
1987) or negative (Dec 1988). This is the effect of the non-linear effects 
that are parametrized by the A^'^ term in Eqs (|I],^Bp. In cold conditions, 
the thermocline has disappeared in our model, so that a Kelvin wave in the 
second baroclinic mode, varying the thermocline depth, is not very effective 
at changing SST. During a warm epsiode the thermocline is so deep that the 
effect of such a wave is also smaller. The effect of a Kelvin wave in the first 
baroclinic mode is not sensitive to the thermocline depth. The effect of these 
nonlinearities in the positive feedback is larger for the NIN03 index than for 
the NIN03.4 index, as the thermocline is shallower in the NIN03 region and 
hence variations of its depth have a larger effect. 

The approximation in the delayed oscillator equation that the Bjerknes 
feedback is instantaneous does not seem to be a very good one, especially 
for the NIN03 index. Part of the delay is due to the artificial 1-month delay 
of the statistical atmosphere, but another 2-3 months are due to the time it 
takes for a Kelvin wave generated at the date line to reach the signal region; 
before this time local negative effects play a role. This effect is also smaller 
for the NIN03.4 index. 

The negative feedback consists of an overlapping train of n = 1 Rossby 
waves that are reflected off different parts of the coast, instead of an idealized 
peak. Combined with the lower effective group velocity shown in Fig. ^ we 
obtain negative influences from 6 to about 18 months prior to the measure- 
ment. Higher order Rossby waves do not contribute at these timescales. 

Due to the normalization uncertainties it is difficult to gauge the extend 
to which the summed positive feedback exceeds one, and compare this quan- 
titatively with the negative feedback terms, i.e., to determine the parameters 
in Eq. (|l]) or Eq. ( p6D in our model. However, it is clear from the figure that 
the negative terms do not depend on the initial state in the same non-linear 
way as the positive terms. The negative feedback is much stronger than ex- 
pected for the December 1988 experiments. We suspect that this is due to 
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the non-linear interactions between the different barochnic modes that we 



observed in section |4.2| . The weak positive response was due to the small 
amplitude at the surface of the first baroclinic mode adjoint Kelvin wave 
that is excited. On the other hand, the adjoint Rossby waves that determine 
the negative feedback propagate in the usual second baroclinic mode in the 
thermocline, which has a larger amplitude at the surface for the same energy 
content and is thus more easily excited by a wind anomaly. At the reflec- 
tion in shallower waters at the western boundary some energy seems to be 
transferred down from Rossby waves in the second baroclinic mode to Kelvin 
waves in the deeper first mode. 



6 Summary and conclusions 

We studied the causes of variations of the sea surface temperature in the 
eastern equatorial Pacific with an adjoint OGCM and statistical atmosphere. 
The adjoint gives us first-order estimates of the sensitivity of the NIN03 
and NIN03.4 indices of this temperature to forcing fields and state fields 
some time earlier. We generated sensitivity fields to these indices at the end 
of October 1987, December 1987 and December 1988; respectively warm, 
transition and cold phases of ENSO. In these experiments we used prescribed 
fluxes to drive the ocean model HOPE. 

The sensitivities either damp out exponentially (heat flux, previous SST, 
zonal surface current) with a time scale of about one month, or propagate as 
equatorial Kelvin and Rossby waves (wind stress, sea level). Of course, with- 
out a full atmosphere model we cannot observe small-scale coupled modes. 

The temperature of the NIN03 region can be affected by a Rossby wave 
coming from the eastern coast. In our model some of this influence is re- 
flected back — an equatorial Kelvin wave impinging on the eastern coast can 
generate such a Rossby wave, albeit not very efficiently. Another way the 
temperature can be changed is by the arrival of a Kelvin wave from the west; 
these are clearly visible in the sensitivities to the wind stress and sea level. 
This Kelvin wave can be generated by Rossby waves reflecting off the coast 
of New Guinea and the Philippines, either directly or with a few months' 
delay at the western boundary area. These delays, combined with the east- 
ern reflections, give rise to an effective group velocity of the sensitivities that 
is lower than the speed of the individual Kelvin and Rossby waves. Higher 
order adjoint Rossby waves are observed, but do not have a strong enough 
wind stress sensitivity field to be able to influence the eastern Pacific. These 
waves sum to a sensitivity to the thermocline depth in the warm pool off 
the equator, but there is no sensitivity to the zonally-averaged thermocline 
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depth along the equator at long timescales. 

The speed of the adjoint Kelvin waves depends very much on the state 
of the ocean in our model at the measurement region. Note that we do 
not simulate Kelvin waves as they occur in the ocean, but as they influence 
the NIN03 index. Normally these propagate in the second baroclinic mode 
(with a zero at the thermocline), with c ~ 2ms~^ in the western and central 
Pacific. However, in the cold tongue during the cold phase (La Nina), the 
thermocline is so shallow that the first baroclinic mode is excited, with a 
speed c ~ 3 ms~^ in the central Pacific. These different modes mix at the 
reflection in shallow waters at the west coast; the reflected adjoint Rossby 
waves have speeds compatible with the second baroclinic modes (the third 
for higher order Rossby waves). These results should be compared to an 
observed speed of 2.4±0.3ms~^, independent of longitude, from the TOGA- 
TAO array, which is ascribed to a mix of the first and second modes. 

We used an adjoint statistical atmosphere model to convert the sensitivi- 
ties into fluxes (mainly zonal wind stress) into sensitivities to the SST. Next 
we approximated these fields by sensitivities to a scalar index: the sensitivity 
to the index region. This underestimates the size of the influence by a factor 
of roughly four, but does not change the shape as a function of time very 
much, due to the high correlations of the neglected terms with the included 
ones. 



The resulting response function (Fig. 13) shows the essence of the delayed 
oscillator mechanism. Going back from the measurement time, the first few 
months exhibit the positive feedback of the change in trade winds reinforcing 
the SST anomaly. At about half a year before the target date the feedback 
changes sign as Rossby waves reflect off the western coast as Kelvin waves 
with opposite polarity. Both these feedback loops are spread out in time. 
Kelvin waves are excited most easily in the central Pacific and take about 
two months to travel to the signal region, and the reflections at the western 
side occur at different longitudes so that the influence of Rossby waves is 
spread out over half a year as well. 

We also observe that the positive feedback is much stronger when the 
system is in transition (Dec 1987), than when it is in a warm (Oct 1987) or 
cold (Dec 1988) phase, pointing to a role of non-linear terms in restraining 
the delayed oscillator. We do not observe this effect in the delayed response. 

In this study we have shown that the adjoint HOPE model gives a good 
approximation of the sensitivity of the full non-linear model more than a 
year back in time. To our surprise most features of the sensitivity fields 
in this primitive equation OGCM can be described by linear long-wave dy- 
namics, modified by the background state. Zonal advection and equatorial 
thermocline recharge mechanisms are not visible. The western-boundary re- 
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flections give rise to a feedback structure similar to tfie one postulated in tlie 
delayed-oscillator meclianism. Unfortunately we cannot measure whether the 
strength of the feedbacks is enough to sustain oscillations, i.e., whether other 
mechanisms, such as stochastic forcing, are needed to excite ENSO episodes. 

Singular vector analysis of simpler but coupled models can be compared 
to the results presented here if the evolved perturbation pattern is similar to 
an ENSO pattern and the difference between our uncoupled response func- 
tion and a coupled one are (qualitatively) accounted for (cf. the discussion 
of Eq. (|25|) ) . The 3- month singular vectors of [Moore and Kleeman (1997a| Jb|) 
are dominated by dynamical processes in the atmosphere, and do not cor- 
respond to our oceanic sensitivities, although their one-month results could 
be compared to our sensitivity to wind fields at two months' lead time (the 
difference is due to the more westerly location of their signal pattern); both 
indicate a high sensitivity to westerly wind bursts. In pCue et al. (1997a| ) the 
final patterns have a typical ENSO signature and should be comparable to 
our sensitivities at 6 months' lead time. They find a predominance of lo- 
cal forcing due to a reduction of upwelling. This mode seems to correspond 
to an iteration of our local forcing in Fig. |l2|a. A new norm (pCue et al.J 
1997b| ) adds to this a wind field comparable with our remote-forcing result 



(Fig. |TO|d), which is dominated by variability of the thermocline depth. This 
signal includes a reflection at the western boundary. 
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A Technical details of the adjoint model 

Most of the adjoint code was generated by automatic differentiation of the 
Fortran code of the HOPE and statistical atmosphere models by the (T)AMC 
( [Giering, 1994| , |1997| ). The exceptions were the horizontal and vertical diffu- 



sion and viscosity computations in HOPE. These involve the inversion of 
a tridiagonal matrix, the adjoint can be written much simpler than the 
automatically-generated algorithm. We have not checked the code for vec- 
torization as all experiments reported here were performed on a scalar CPU 
(MIPS RIOOOO). 

The adjoint model depends both on the adjoint fields d/dx and the 
fields of the reference trajectory x for non-linear functions Ai. The latter 
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fields are saved from tlie forward run of tfie model to be used in the backward 
run of the adjoint model in a three-layer scheme. To compute the adjoint 
of one day (12 ocean time steps), the model is run forward over this period 
while storing fields needed in non-linear transformations on a direct-access 
file. Saving these in single precision, and eliminating the multiple storage of 
the slowly-varying temperature and salinity fields per time step one obtains 
a 250 MB storage file (21 3D fields of 2*97*61*20 numbers per time step 
times 12 steps per day). These daily runs are restarted from restart files 
generated by a one-month run. Each (modified) restart file takes 8 MB (8 
3D fields). The monthly restart files are generated by a top-level run. In all 
the storage requirements are about 500 MB disk space. The run-time size of 
the program is 130 MB. 

Apart from end effects (the last month can be skipped on the highest 
level) the running time for an adjoint run is three forward runs (total, month 
and day) and one backward run. The latter includes many forward recom- 
putations, and takes about 2-3 times longer than a forward run. The total 
cost of a full adjoint run is therefore 5-6 times the cost of a forward run in 
CPU time. To this one should add the cost of writing and reading 250 MB 
temporary storage per simulated day. 

We verified whether finite perturbations of the initial state fields did in 
fact give results in agreement with the adjoint model, which is only strictly 
valid for infinitesimal perturbations. For small random perturbations (1-6 
wavelengths in x and y and 0-1 in z tapering off with depth) of all initial 
state variables except temperature and salinity we obtain good agreement 
with the derivative up to more than 6 months. We did not include variations 
of the temperature and salinity fields as these give rise to different convection 
patterns, mainly in the mixed layer at high latitudes. The derivative to 
the most important input variable, Tx, is validated more systematically in 
Appendix |B[ 



B Validation 

Although the linear dynamics we found in section | looks plausible, it remains 
to be verified that non-linear effects do not in fact dominate in the model 
(let alone in nature). For this purpose we compared the influence of zonal 
wind fields on the NIN03 index at the end of 1987 with the predictions of 
the adjoint model, cf. Eq. (0) 

5Ns ^ ^6tx . (28) 
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Figure 14: The quality of the adjoint model 1986-1987 run compared with 
perturbed forward runs. On the bottom plane a space-time diagram of the 
sensitivity to the zonal wind field, dN^/dTx, along the equator is sketched 
(dashed contours). The groups of four bars represent the quality of the ad- 
joint field, checked with perturbed wind stresses of St^ = —0.05, —0.01, 
-1-0.01, -1-0.05 Nm~^ respectively. The position of the bars indicates the lon- 
gitude and date of the perturbation. The height gives the quality, defined as 
q = 1 — 2/tt a,Tcta,n\l — 6N^/{6TxdN3/dTx)\. It is 1 for perfect agreement, 0.6 
when the ratio is 0.5 or 1.5, and it tends to zero for very bad agreement 



All deviations are relative to a forward run with prescribed fluxesQ. For the 
extra wind field St^ we chose a longitude with a large influence, the size 
was two grid points at latitude ±1.4° with a zonal extent of 20°; the zonal 
edges are smoothed over ±5°. The duration of the storm was matched to 
the sampling of the derivative fields, 8 days, and its intensity was ±0.01 and 
±0.05 Nm~^, resulting in predicted deviations of the NIN03 index of roughly 
0.05 K (1 month) to 0.01 K (1 year) for the higher wind stress. 

In Fig. |l^ we present the results of these experiments. In the initial ad- 
joint Kelvin wave (Sep-Dec) the ratios almost always are within the range 
0.8-1.3. Smaller wind fields give better results, as the error is dominated 

'^Due to an oversight, relaxed to climatology rather than observed SST; this should not 
influence the findings. 
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by non-linear effects. The exceptions are in September, when the reffected 
Rossby waves inffuence the picture. The previous two months (Jul-Aug), 
when the reffections are building up and the perturbations are near the west- 
ern coast, the adjoint performs poorly. For higher wind fields the ration lies 
between 0.2 and 2. These perturbations send waves of both signs across the 
Pacific, so the final perturbation of the NIN03 index is a sum of positive and 
negative patches that is very sensitive to details. 

The situation improves when the adjoint Rossby waves are clearly estab- 
lished (Jan — Jun). When the signal is not too weak the prediction for the 
stronger winds is accurate to within a factor 0.4-1.6 (it tends to be a bit 
too low). The error is dominated by noise, so that the lower wind fields give 
worse results. In the previous year we only investigated the stronger wind 
fields. One sees the same pattern: the well-defined first Rossby wave gives 
reasonable agreement (although the derivative again underestimates the real 
effect by a factor 1.5 on average). However, the effect of a wind field outside 
the strongest wave is predicted poorly by the adjoint. 



Date 




longitude 




5NIN03 


adjoint 


ratio 








[Nm-2] 


[K] 


[K] 




Oct 19^ 


17 


160°-180° 


0.2 


0.295 


0.199 


1.48 








-0.2 


-0.173 


-0.199 


0.86 


Jun 19^ 


17 


160°-180° 


0.2 


-0.041 


-0.069 


0.60 








-0.2 


0.062 


0.069 


0.90 


Jan 19^ 


17 


225°-245° 


0.2 


-0.048 


-0.048 


1.00 








-0.2 


0.060 


0.048 


1.25 


Aug 19^ 




160°-180° 


0.2 


0.000 


-0.019 


-0.02 








-0.2 


0.012 


0.019 


0.64 



Table 1: The effect of 8-day stronger anomalous wind bursts over 20° x 4° 
along the equator. 

We also studied the effect of stronger anomalous wind stress fields, 5tx = 
0.2 Nm~^, which is the typical strength of a westerly wind burst ( |Giese and 



Harrison, 1990|) . The results are presented in Table ^ Non- linear effects 
worsen the agreement at short lead times compared to more gentle pertur- 
bations, but at longer lead times the results are comparable. 

In all, it seems we can trust the sensitivities of the adjoint model most 
of the time to within a factor 2, even to one year back. This does not hold 
when the expected deviation is small, or when the perturbation is near the 
western coasts so that reflected and direct waves interfere. 
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